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Abstract. - It is shown that non-Markovian master equations for an open system which are local 
in time can be unravelled through a piecewise deterministic quantum jump process in its Hilbert 
space. We derive a stochastic Schrodinger equation that reveals how non-Markovian effects are 
manifested in statistical correlations between different realizations of the process. Moreover, we 
demonstrate that possible violations of the positivity of approximate master equations are closely 
connected to singularities of the stochastic Schrodinger equation, which could lead to important 
insights into the structural characterization of positive non-Markovian equations of motion. 
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Relaxation and decoherence phenomena in open quan- 
tum systems [1] can often be modelled with sufficient ac- 
curacy by a quantum Markov processes in which the open 
system's density matrix is governed by a relatively sim- 
ple quantum Markovian master equation with Lindblad 
structure [2,3]. However, non-Markovian quantum sys- 
tems featuring strong memory effects play an increasingly 
important role in many fields of physics such as quantum 
optics [4] , solid state physics [5] , and quantum information 
science [6]. Further applications include non-Markovian 
extensions of quantum process tomography, quantum con- 
trol [7], and quantum transport [8]. 

The non-Markovian quantum dynamics of open systems 
is characterized by pronounced memory effects, finite re- 
vival times and non-exponential behavior of damping and 
decoherence, resulting from long-range correlation func- 
tions and from the dynamical relevance of large correla- 
tions and entanglement in the initial state. As a conse- 
quence the theoretical treatment of non-Markovian quan- 
tum dynamics is generally extremely demanding, both 
from the analytical and from the computational point of 
view [9] . Even if one is able to derive an appropriate non- 
Markovian master equation or some other mathematical 
formulation of the dynamics, the numerical simulation of 
such processes turns out to be a very difficult and time- 
consuming task, especially for high-dimensional Hilbert 
spaces. 

From classical physics it is known that Monte Carlo 
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techniques provide efficient tools for the numerical simu- 
lation of complex systems. This fact was one of the moti- 
vations to introduce the Monte Carlo wave function tech- 
nique [10-12] which provides efficient quantum simulation 
techniques in the regime of Markovian dynamics. Sev- 
eral generalizations of the Monte Carlo approach to non- 
Markovian dynamics have been developed which are based 
on suitable extensions of the underlying reduced system's 
Hilbert space [13-16]. 

Recently, an efficient alternative simulation algorithm 
for the treatment of non-Markovian open system dynam- 
ics has been proposed [17] that does not require any exten- 
sion of the state space. The purpose of the present paper 
is to develop a mathematical formulation of this algorithm 
in terms of a stochastic Schrodinger equation (SSE) in the 
open system's Hilbert space. We demonstrate that this 
formulation gives rise to a new type of piecewise deter- 
ministic quantum jumps process (PDP). 

Quantum master equations are often derived from an 
underlying microscopic theory by employing some ap- 
proximation scheme. An appropriate scheme is the 
time-convolutionlcss (TCL) projection operator technique 
which leads to a time-local first-order differential equation 
for the density matrix [18-20]. ft will be shown that TCL 
master equations allow a stochastic unravelling of the form 
developed here. Generally, the use of a certain approxi- 
mation technique may lead to violations of the positivity 
of the master equation. We demonstrate that positivity 
violations are closely linked to singularities of the SSE at 
which the stochastic process breaks down. Hence, a great 
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advantage of the present stochastic formulation consists 
in the fact that it naturally prevents the generation of 
unphysical solutions and that it could thus lead to impor- 
tant insights into the structural characterization of posi- 
tive non-Markovian evolution equations. 

The most general structure of the TCL master equation 
is given by 1 



d , , 



Ctp(t) = -i[H{t),p{t))+V t p{t), (1) 



where 



T> t p = J2 A m(t) 



c m (t) P cl(t)-±{c'L(t)c m (t),p} 



The time-dependent generator Ct consists of a commu- 
tator term describing the unitary part of the evolution 
and a dissipator T> t . The latter involves a summation over 
the various decay channels labelled by m with correspond- 
ing time-dependent decay rates A m (t) and arbitrary time- 
dependent system operators C m (t). 

In the simplest case the rates A m as well as the Hamil- 
tonian H and the operators C m are assumed to be time- 
independent. Equation (1) then represents a master equa- 
tion in Lindblad form [2,3] which generates a semigroup 
of completely positive dynamical maps known as quan- 
tum Markov process. For arbitrary time-dependent oper- 
ators H{t) and C m (t), and for A m (t) > the generator 
of the master equation (1) is still in Lindblad form for 
each fixed time t and leads to a two-parameter family of 
completely positive dynamical transformations [21] which 
may be referred to as time-dependent quantum Markov 
process [22] . An entirely different situation emerges if one 
or several of the A m (t) become temporarily negative which 
expresses the presence of strong memory effects in the re- 
duced system dynamics. The process is then said to be 
non-Markovian. Of course, the physical interpretation of 
the master equation requires that it preserves the positiv- 
ity of the density matrix p. The formulation of general 
mathematical and physical conditions that guarantee the 
preservation of positivity is, however, an unsolved problem 
of central importance in the field of non-Markovian quan- 
tum dynamics. We emphasize that the emergence of tem- 
porarily negative A rn (t) in the master equation is a natural 
phenomenon in the non-Markovian regime which does in 
general not imply that the complete positivity of the corre- 
sponding quantum dynamical map is violated. An exam- 
ple is discussed in Ref. [1] where the exact non-Markovian 
master equation of an analytically solvable model is con- 
structed. 

The fundamental difference between Markovian and 
time-dependent Markovian processes on the one hand 
and non-Markovian processes on the other hand can also 
be seen very clearly if one attempts to apply the stan- 
dard stochastic formulations to the master equation (1). 



1 This statement is a direct consequence of Lemma 2.3 of Ref. [2]; 
for various alternative forms see [1,16]. 



For both a Markovian and a time-dependent Markovian 
dynamics the standard unravelling through a stochastic 
quantum jump process can indeed be applied. This means 
that in both cases one can formulate an appropriate PDP 
for the reduced system's state vector \ip(t)) in such a way 
that the expectation value 



p(t) = E[\m)(m\ 



di>P[\^),t]\^)(iP\ (2) 



satisfies the master equation (1). Here, we have expressed 
the expectation value E in terms of an integration over the 
Hilbert space of states of the open quantum system with 
the unitarily invariant volume clement dip = DiJjDip* , and 
introduced the corresponding probability density func- 
tional P[|-!/)),£] which is defined as the probability density 
of finding at time t the state vector \ip) [1]. However, the 
essential feature of a non-Markovian dynamics is the tem- 
porary appearance of negative decay rates. The use of 
the standard unravellings unavoidably leads in this case 
to negative jump probabilities, which clearly indicates the 
decisive difference between Markovian and non-Markovian 
quantum processes. 

To account for the sign of the decay rates we decom- 
pose A m (t) into a positive and a negative part defined by 
A±(() = |[|A m (i)| ± A m (t)]. The master equation (1) 
can then be written in the form 



dt P 



-i[H(t),p] (3) 

X>+(t) \c k (t) P ct(t) - \ {cl(t)c k (t) lP } 

k L 

X>rw \ci(t) P ci{t) - \ \c\{t)c x {t\ P } 



Iii order to better distinguish the positive and the nega- 
tive channels we label the former by an index k and the 
latter by an index I. Note that A+(<) > and Af(t) > 
and that for Markovian or time-dependent Markovian pro- 
cesses we have A ( ~ (t) = 0. 

We can now formulate the central result of this paper. 
Namely, the stochastic Schrodinger equation given by 
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-iG(t)\^(t))dt 

? [\mmm\ ~ m)) \ dNk {t) 
J2 I dip'[W)-m))]dN^,{t) (4) 



yields an unravelling of the master equation (3) through a 
non-Markovian piecewise deterministic process. The first 
term on the r.h.s. of Eq. (4) represents the normalized 
deterministic drift of the process which is generated by 
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The instantaneous and random quantum jumps are de- 
scribed by the second and the third line of Eq. (4). The 
quantities dN^(t) and dN^,(t) are random Poisson in- 
crements satisfying the relations 



dN+(t)dN+(t) 
dN^,(t)dN^„(t) 
dN+(t)dN^,,(t) 



6 kk ,dN+(t), 

5 w t{W)-W'))dN^,{t), 
0, (5) 



and having expectation values 

E[dN+(t)] = A+(t)(m\Cl(t)C k (t)\m)dt, (6) 

E[dNr,,(t)} = AY(t) ^Q4 (ij'\cjm(tM) 
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\m)w)\\ 



dt. (7) 



Here, the delta functional on Hilbert space is defined by 
~ \ipo))F[\i>}} = F[\1>o)], where F[\1>)] is an ar- 
bitrary smooth functional. 

The physical meaning of the properties in (5) is that 
there cannot be two or more jumps simultaneously in a 
given realization of the process and in a given moment of 
time. Suppose first that the dynamics is Markovian or 
time-dependent Markovian. We then have dN^,(t) = 
and the SSE (4) reduces to the stochastic differential equa- 
tion of the standard PDP unravelling. According to the 
second line of Eq. (4) the quantum jumps are represented 
by an instantaneous change of the state vector, 



c h (t)W)) 



\\Ck(t)Mt))W 

and by virtue of Eq. (6) this jump occurs at the rate 



A+(t)m)\ct(t)c k (t)m))- 



(8) 



The term in the third line of the SSE (4) describes the 
negative channels which arc crucial for the unravelling of 
non-Markovian dynamics. The corresponding jumps are 
given by instantaneous transitions from the actual state 
\ip(t)) to some state \ip')- To account for all possible tar- 
get states of the negative channel jumps we perform in this 
term an integration over \ip'). According to the delta func- 
tional in Eq. (7) the target state \ip') of the possible jumps 
is related to the source state \ip) by \tp) = C;|i/ , '}/||C;|V , ')ll- 
Hence, negative jumps correspond to a reversal of cer- 
tain positive jumps, obtained by interchanging the role 
of source and target state. The quantity dN^,(t) is the 
Poisson increment for the negative jumps via channel I. 
From Eq. (7) we infer that the state vector \ip) can per- 
form a jump to a state vector in some volume element dip' 
of Hilbert space around \ip') with the rate (for simplicity 
we omit the time arguments) 



\\CiW)\\ 



dip. 



In an ensemble of realizations of the process the quantity 
P[\ip')]dip' / P[\ip)]dip can be interpreted as N'/N, where 
N' is the number of realizations in volume element dip' 
and N is the number of realizations in clement dip. Then 
we can identify 5 (\ip) - Q\ip'} /\\Ci\ip')\\) dip = 1. Hence, 
the negative channel jumps from \ip) to \ip') occur at the 
rate 

T-=Ai^(i ) '\C}C l W). (9) 

The comparison with the rate (8) for positive jumps shows 
two crucial differences. First, the rates for the positive 
jumps is proportional to the expectation value of Cj,C k in 
the source state, while the rates for the negative jumps 
is proportional to the expectation value of CjCi in the 
target state. Hence, again the role of source and target 
state have been interchanged. Second, the negative jump 
rates carry an additional factor of N'/N, the ratio of the 
number of ensemble members in the target state to the 
number of members in the source state. Note that due to 
the presence of this factor the SSE (4) is not a stochastic 
differential equation in the usual sense because the expec- 
tation values of the random increments (7) depend explic- 
itly on the full probability density. To determine these 
increments at a certain time t one has to know the proba- 
bility density P[|^) , t] . Within a numerical simulation this 
is achieved by propagating simultaneously an ensemble of 
realization from which t] can then be estimated sclf- 

consistently. As an important consequence and as a result 
of the non-Markovian character of the dynamics we thus 
find certain correlations between different realizations of 
the process. 

It may seem at first sight that the correlations between 
the realizations require that a huge number of realizations 
of the process has to be generated simultaneously in order 
to obtain a sufficiently accurate estimate for the probabil- 
ity density. However, when the realizations of the process 
are generated on a computer there is no need to have iVj 
identical copies of the state \ipi) to obtain P \\ipi)]. It is suf- 
ficient to have only a single copy of \ipi) and to keep track 
of the corresponding integer number JV». This allows to 
optimize the numerical implementation of the process and 
to perform simulations in an efficient way [17]. 

To prove that the expectation value (2) for the pro- 
cess obtained from the SSE (4) indeed satisfies the master 
equation (3) we start from 

d(ivwi) = + w(dip\ + midipi. 

Taking the expectation value of this relation, expressing 
the increments \dip) through the SSE (4), and using the 
properties (5) we find 

dp = - l \H,p}dt-J2^ 1 {CiC m ,p}dt 

m 

+ J2 A rnE[(lP\ClC m \lP)\lP)(lP\]dt 
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Using here the expectation values of the increments from 
Eq. (7) one immediately obtains the required master 
equation (3). Hence, we have proven the validity of 
SSE (4) which is the central result of the paper. It 
is worth mentioning that while there exists stochastic 
Schrodingcr equations of diffusion type for non-Markovian 
systems [23,24], to the best of our knowledge our SSE 
is the first representation through a stochastic quantum 
jump process in the reduced system's Hilbert space. 

It is important to note that the expectation value for 
dN^,(t) in Eq. (7) is not well defined when the denom- 
inator becomes equal to zero, i. e. P[|-0),i] = 0, where 
!</,) = Ci\ip')/\\Ci\ip')\\, or alternatively N = in Eq. (9). 
The stochastic process breaks down at this point since 
there exists an open negative channel but there are no re- 
alizations which are in the source state of the correspond- 
ing jump. The formulation of general conditions on the 
structure of the master equation (1) that ensure the ab- 
sence of such singularities in the corresponding SSE (4) 
seems to be a difficult problem. However, it is quite easy 
to demonstrate that a breakdown of the process necessar- 
ily takes place if the master equation violates positivity. 
In fact, within the stochastic formulation developed here 
the density matrix p(t) of the open system is given by 
the expectation value (2) which represents, by the very 
construction, a positive matrix. Therefore, if the mas- 
ter equation violates positivity at some point of time the 
stochastic dynamics must necessarily cease to exist. The 
present method thus signals the point of violation of pos- 
itivity of the density matrix. 

To prove this statement we denote the state space, 
i. e., the set of all density matrices of the open system 
by S. Let us assume that the master equation (1) vi- 
olates positivity. Hence, there is an initial state p(Q) 
and a corresponding solution p{t) of the master equation 
which leaves the state space S after some point of time 
t = to. At this point p(to) = po reaches the boundary 
of S. Let A(t) = (ip(t)\p(t)\ip(t)) be the lowest eigenvalue 
of p(t) with corresponding eigenvector \ip(t)). Lying on 
the boundary, po must have at least one zero eigenvalue 
with corresponding eigenvector |cpo) = |<p(^o))j i- e., we 
have A(£o) = ( l Po\po\ l Po) = 0. Hence, an appropriate con- 
dition implying the violation of positivity is given by the 
inequality A (to) < 2 . The Hellman-Feynman theorem 
yields 

A(*o) = (<A)|p(io)|<A)} = (^olAoPol^o), 
and we find the following condition for the violation of 



(</?o|£t Po|<A)} < 0. 



(10) 



Consider now an ensemble representation of po that is 
generated through the SSE (4): po = '}2, i Pi\4'i){4'i\ with 
(ipi\ipi) = 1, Pi > and YliPi = 1- We- then have 

(tpo\po\tpo) = ^p 4 |(^o|^:)| 2 = 0. 



It follows that \(po) is orthogonal to all members of the 
ensemble, i. e., (<po|'0i) = 0. Evaluating condition (10) 
one therefore finds 



(<po\£t o po\tp ) 



p i Am(t )\(ipo\C m \'ip l }\ 2 < 0. 



2 It is assumed here that A(to) does not vanish, which is obviously 
the generic case. 



Hence, there must exist indices m and i such that 
A m (to) < and (</2o|C m |V>i) ^ 0. It follows that C m \tpi) 
has a nonzero component in the direction of \tpo) and that 
the state vector \ijj) = C m \ipi) /\\C m \i}i)\\ does not belong 
to the ensemble In other words, P to] = 0. We 

conclude that the point of violation of positivity implies 
the breakdown of the SSE (4) because there exists an open 
channel with negative rate while the probability of being 
in the source state of the corresponding jump vanishes. 

Of course, the formal mathematical solution of the mas- 
ter equation (1) does not halt at the point of time when 
the positivity is lost: The dynamics continues to reduce 
occupation probability of a given state beyond the zero- 
point. However, the evolution given by the SSE (4) stops 
at the zero-point since the number of realizations in a 
given state cannot, by construction, have negative val- 
ues. Thus, the stochastic process developed here identi- 
fies the point of time where the master equation loses the 
positivity, preventing excursions to unphysical solutions. 
While the stochastic unravelling of the master equation is 
in general not unique, we expect that the connection be- 
tween a breakdown of the positivity and a singularity of 
the SSE holds for all stochastic representation of the form 
constructed here. 

In conclusion, we have derived a piecewise deterministic 
process which describes the dynamics of non-Markovian 
systems. The stochastic Schrodingcr equation constructed 
reveals the fundamental mathematical and physical differ- 
ence between time-local master equations which appear 
with positive and with negative rates. The corresponding 
Poisson increments have a distinct structure and the nega- 
tive rate process clearly shows how non-Markovian effects 
are manifested. 

Markovian and non-Markovian processes are widely 
used for the modelling of dynamical systems in many ar- 
eas of physics, chemistry and biophysics. Our results indi- 
cate how to treat master equations with negative rates and 
memory effects also for classical systems. In fact, the stan- 
dard simulation algorithm for a classical Markovian mas- 
ter equation corresponding to the stochastic wave function 
method is known as Gillespie algorithm [25]. The method 
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proposed here could therefore lead to the development of [24] Diosi L., Gisin N. and Strunz W. T., Phys. Rev. A, 58 
an efficient non-Markovian generalization of the Gillespie (1998) 1699. 

algorithm and thus opens the way to many further studies [25] Gillespie D. T., J. Phys. Chem., 81 (1977) 2340. 
in the dynamics of complex system. 
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